# -*- coding: utf-8 -*-
"""
Created on Tue Apr 09 11:54:04 2013

@author: Maxim
"""
import aircraft
import FlightConditions as FC
import aerodynamics as aero
import matplotlib as mpl
import numpy as np

def experimental_data():
    CL = np.array([-0.5905,-0.4236,-0.2512,-0.069,0.1195,0.3064,0.491,0.6706,
                   0.8343,0.9837,1.0843,1.1531,1.217,1.1183,0.9145,0.9078,0.915,
                   0.9298])
    CD = np.array([0.0762,0.0451,0.0346,0.0291,0.0284,0.0311,0.0371,0.0468,0.0596,
                   0.076,0.0969,0.1247,0.1627,0.2557,0.3492,0.4044,0.455,0.5065])
    alpha  = np.array([-10.1532,-8.0954,-6.0376,-3.9635,-1.8992,0.1691,2.2383,4.311,
                       6.3628,8.4257,10.4613,12.4963,14.5141,16.4196,18.3688,
                       20.4194,22.3854,24.3961])
    return CL, CD, alpha

def aero_valid():
    ac = aircraft.load('V0510')
    fc = aero.flightConditions()
    CG = ac.get_CG(False)
    I = ac.get_inertia(False)
    M = ac.get_mass_total(False)
    V = 34.0
    rho = 1.2255
    CD0 = ac.get_drag(V,0.0)
    fc.addTrimmedFlightCondition('trim',M, CG, I, V, rho)
    aero3D = aero.aero3d_VLM(ac)
    rslt = aero3D.runVLM(fc)
    
    alpha = np.linspace(-10.,20.,20)
    CLexp, CDexp, alphaExp = experimental_data()
    
    CLa = rslt.results[0].a * (np.pi/180.0)
    CL0 = rslt.results[0].CL0
    k   = rslt.results[0].k
    rslt.results[0].display()
    CL = CL0 + alpha*CLa
    CD = CD0 + CL*CL*k

    plt = mpl.pyplot
    
    font = {'family' : 'normal',
            'weight' : 'normal',
            'size'   : 20}
    
    mpl.rc('font',**font)    
    
    plt.figure(1)
    plt.hold(True)
    plt.grid(True)
    plt.plot(CDexp,CLexp,'bo-')
    plt.plot(CD,CL,'r-')
    plt.ylabel('lift coefficient')
    plt.xlabel('drag coefficient')
    plt.axis([0,0.3,-1.0,1.5])
    plt.legend(['wind tunnel','calculated'],'lower right')
    
    plt.figure(2)
    plt.hold(True)
    plt.grid(True)
    plt.plot(alphaExp,CDexp,'bo-')
    plt.plot(alpha,CD,'r-')
    plt.xlabel('angle of attack,deg')
    plt.ylabel('drag coefficient')
    plt.legend(['wind tunnel','calculated'],'upper left')
    
    plt.figure(3)
    plt.hold(True)
    plt.grid(True)
    plt.plot(alphaExp,CLexp,'bo-')
    plt.plot(alpha,CL,'r-')
    plt.xlabel('angle of attack,deg')
    plt.ylabel('lift coefficient')
    plt.legend(['wind tunnel','calculated'],'lower right')
    plt.axis([-10,25,-0.5,1.5])
    plt.show()

def get_drag_polar():
    rsltFile = 'drag_polar.txt'
    
    ac = aircraft.load('V0510')
    fc = aero.flightConditions()
    CG = ac.get_CG(False)
    I = ac.get_inertia(False)
    M = ac.get_mass_total(False)
    V = 34.0
    rho = 1.2255
    CD0 = ac.get_drag(V,0.0)
    fc.addTrimmedFlightCondition('trim',M, CG, I, V, rho)
    aero3D = aero.aero3d_VLM(ac)
    rslt = aero3D.runVLM(fc)
    
    alpha = np.linspace(-10.,20.,20)
    CLexp, CDexp, alphaExp = experimental_data()
    
    CLa = rslt.results[0].a * (np.pi/180.0)
    CL0 = rslt.results[0].CL0
    k   = rslt.results[0].k
    CL = CL0 + alpha*CLa
    CD = CD0 + CL*CL*k
    
    fid = open(rsltFile,'wt')
    fid.write('Velocity %.2f m/sec\n'%V)
    fid.write('Density %.2f kg/m3\n'%rho)
    fid.write('CD0 %.4f\n'%CD0)
    fid.write('alpha[deg]\tCL\tCD\n')
    for i,a in enumerate(alpha):
        fid.write('%.1f\t%.4f\t%.4f\n'%(a,CL[i],CD[i]))
    fid.close()

if __name__=="__main__":
    aero_valid()